

/*
Function
cblas_cgbmv
Scales a general band matrix, then multiplies by a vector, then adds a vector (single-precision complex).

https://developer.apple.com/documentation/accelerate/1513332-cblas_cgbmv?language=objc

Parameters
Order
Specifies row-major (C) or column-major (Fortran) data ordering.

TransA
Specifies whether to use matrix A ('N' or 'n'), the transpose of A ('T' or 't'), or the conjugate of A ('C' or 'c').

M
Number of rows in matrix A.

N
Number of columns in matrix A.

KL
Number of subdiagonals in matrix A.

KU
Number of superdiagonals in matrix A.

alpha
Scaling factor to multiply matrix A by.

A
Matrix A.

lda
Leading dimension of array containing matrix A. (Must be at least KL+KU+1.)

X
Vector X.

incX
Stride within X. For example, if incX is 7, every 7th element is used.

beta
Scaling factor to multiply vector Y by.

Y
Vector Y.

incY
Stride within Y. For example, if incY is 7, every 7th element is used.

Discussion
Computes alpha*A*x + beta*y, alpha*A'*x + beta*y, or alpha*conjg(A')*x + beta*y depending on the value of TransA.

*/

use common::gsl_cblas_h::{CBLAS_ORDER, CBLAS_TRANSPOSE};

use common::complex::complex::Complex;
use crate::cblash_h;

pub fn cblas_cgbmv (order: &CBLAS_ORDER , TransA: &CBLAS_TRANSPOSE,
    M: i32, N: i32, KL: i32, KU: i32,
    alpha: &Complex, A: &Complex, lda: i32, X: &Vec<Complex>,
    incX: i32, beta: &Complex, Y: &mut Vec<Complex>, incY: i32) {
    //#include "source_gbmv_c.h"

    // let mut i, j;
    // INDEX lenX, lenY, L, U;
  
    let alpha_real = alpha.real;
    let alpha_imag = alpha.imag;
  
    let beta_real = beta.real;
    let beta_imag = beta.imag;
  
    // CHECK_ARGS14(GBMV,order,TransA,M,N,KL,KU,alpha,A,lda,X,incX,beta,Y,incY);
  
    if M == 0 || N == 0 {
        return;
    }

    
    if alpha.equals_0_0() || beta.equals_1_0() {
        return;
    }
    let (lenX, lenY, L, U) = match TransA {
        CBLAS_TRANSPOSE::CblasNoTrans => {
            (N as usize, M as usize, KL, KU)
        },
        CBLAS_TRANSPOSE::CblasTrans | CBLAS_TRANSPOSE::CblasConjTrans => {
            (M as usize, N as usize, KU, KL)
        }
    };

  
    /* form  y := beta*y */
      let iy = cblash_h::offset(lenY, incY);
    if beta.equals_0_0() {
        // beta_real == 0.0 && beta_imag == 0.0
        for index in iy.indices.iter() {
            let item_y = Y.get_mut(*index).expect("Unreachable");
            item_y.real = 0.0;
            item_y.imag = 0.0;
        }
    } else if !beta.equals_1_0() {
        for index in iy.indices.iter() {
            let item_y = Y.get_mut(*index).expect("Unreachable");
            let y_real = item_y.real;
            let y_imag = item_y.imag;
            item_y.real = y_real * beta_real - y_imag * beta_imag;
            item_y.imag = y_real * beta_imag + y_imag * beta_real;
        }
    }

    let ix = cblash_h::offset(lenX, incX);

    match (order, TransA) {
        (CBLAS_ORDER::CblasRowMajor, CBLAS_TRANSPOSE::CblasNoTrans) | (CBLAS_ORDER::CblasColMajor, CBLAS_TRANSPOSE::CblasTrans) => {
            /*
            if ((order == CblasRowMajor && TransA == CblasNoTrans)
            || (order == CblasColMajor && TransA == CblasTrans))
             */
            /* form  y := alpha*A*x + y */
            //     for index in iy.indices.iter() {
            //         let item_y = Y.get_mut(*index).expect("Unreachable");
            //     let dotR = 0.0;
            //     let dotI = 0.0;
            //     let j_min = if i > L { i - L } else { 0 };
            //     let j_max = lenX.min(i + U + 1);
            //     INDEX ix = OFFSET(lenX, incX) + j_min * incX;
            //     for (j = j_min; j < j_max; j++) {
            //     const BASE x_real = CONST_REAL(X, ix);
            //     const BASE x_imag = CONST_IMAG(X, ix);
            //     const BASE A_real = CONST_REAL(A, lda * i + (L + j - i));
            //     const BASE A_imag = CONST_IMAG(A, lda * i + (L + j - i));
        
            //     dotR += A_real * x_real - A_imag * x_imag;
            //     dotI += A_real * x_imag + A_imag * x_real;
            //     ix += incX;
            //     }
        
            //     REAL(Y, iy) += alpha_real * dotR - alpha_imag * dotI;
            //     IMAG(Y, iy) += alpha_real * dotI + alpha_imag * dotR;
            //     iy += incY;
            // }
        },
        (CBLAS_ORDER::CblasRowMajor, CBLAS_TRANSPOSE::CblasTrans) | (CBLAS_ORDER::CblasColMajor, CBLAS_TRANSPOSE::CblasNoTrans) => {
            /*
            else if ((order == CblasRowMajor && TransA == CblasTrans)
               || (order == CblasColMajor && TransA == CblasNoTrans)) {
             */
            
            // /* form  y := alpha*A'*x + y */
            // INDEX ix = OFFSET(lenX, incX);
            // for (j = 0; j < lenX; j++) {
            //     const BASE x_real = CONST_REAL(X, ix);
            //     const BASE x_imag = CONST_IMAG(X, ix);
            //     BASE tmpR = alpha_real * x_real - alpha_imag * x_imag;
            //     BASE tmpI = alpha_real * x_imag + alpha_imag * x_real;
            //     if (!(tmpR == 0.0 && tmpI == 0.0)) {
            //     const INDEX i_min = (j > U ? j - U : 0);
            //     const INDEX i_max = GSL_MIN(lenY, j + L + 1);
            //     INDEX iy = OFFSET(lenY, incY) + i_min * incY;
            //     for (i = i_min; i < i_max; i++) {
            //         const BASE A_real = CONST_REAL(A, lda * j + (U + i - j));
            //         const BASE A_imag = CONST_IMAG(A, lda * j + (U + i - j));
            //         REAL(Y, iy) += A_real * tmpR - A_imag * tmpI;
            //         IMAG(Y, iy) += A_real * tmpI + A_imag * tmpR;
            //         iy += incY;
            //     }
            //     }
            //     ix += incX;
            // }
        },
        (CBLAS_ORDER::CblasRowMajor, CBLAS_TRANSPOSE::CblasConjTrans) => {
            /*
            } else if (order == CblasRowMajor && TransA == CblasConjTrans) {
             */
            
            // /* form  y := alpha*A^H*x + y */
            // INDEX ix = OFFSET(lenX, incX);
            // for (j = 0; j < lenX; j++) {
            //     const BASE x_real = CONST_REAL(X, ix);
            //     const BASE x_imag = CONST_IMAG(X, ix);
            //     BASE tmpR = alpha_real * x_real - alpha_imag * x_imag;
            //     BASE tmpI = alpha_real * x_imag + alpha_imag * x_real;
            //     if (!(tmpR == 0.0 && tmpI == 0.0)) {
            //     const INDEX i_min = (j > U ? j - U : 0);
            //     const INDEX i_max = GSL_MIN(lenY, j + L + 1);
            //     INDEX iy = OFFSET(lenY, incY) + i_min * incY;
            //     for (i = i_min; i < i_max; i++) {
            //         const BASE A_real = CONST_REAL(A, lda * j + (U + i - j));
            //         const BASE A_imag = CONST_IMAG(A, lda * j + (U + i - j));
            //         REAL(Y, iy) += A_real * tmpR - (-A_imag) * tmpI;
            //         IMAG(Y, iy) += A_real * tmpI + (-A_imag) * tmpR;
            //         iy += incY;
            //     }
            //     }
            //     ix += incX;
            // }
        },
        (CBLAS_ORDER::CblasRowMajor, CBLAS_TRANSPOSE::CblasConjTrans) => {
            /*
            } else if (order == CblasColMajor && TransA == CblasConjTrans) {
             */
            
            // /* form  y := alpha*A^H*x + y */
            // INDEX iy = OFFSET(lenY, incY);
            // for (i = 0; i < lenY; i++) {
            //     BASE dotR = 0.0;
            //     BASE dotI = 0.0;
            //     const INDEX j_min = (i > L ? i - L : 0);
            //     const INDEX j_max = GSL_MIN(lenX, i + U + 1);
            //     INDEX ix = OFFSET(lenX, incX) + j_min * incX;
            //     for (j = j_min; j < j_max; j++) {
            //     const BASE x_real = CONST_REAL(X, ix);
            //     const BASE x_imag = CONST_IMAG(X, ix);
            //     const BASE A_real = CONST_REAL(A, lda * i + (L + j - i));
            //     const BASE A_imag = CONST_IMAG(A, lda * i + (L + j - i));
        
            //     dotR += A_real * x_real - (-A_imag) * x_imag;
            //     dotI += A_real * x_imag + (-A_imag) * x_real;
            //     ix += incX;
            //     }
        
            //     REAL(Y, iy) += alpha_real * dotR - alpha_imag * dotI;
            //     IMAG(Y, iy) += alpha_real * dotI + alpha_imag * dotR;
            //     iy += incY;
            // }
        },
        _ => panic!("BLAS_ERROR: unrecognized operation"),
    }
    
}